*------------------------------------------------------------------------------*
* Dependent variables: a11, b9 and C6a
* Observable implications DVs: D6 (ordered logit), D10 (logit), D11 (ordered logit)
*------------------------------------------------------------------------------*


*------------------
* Prepare working space
*------------------

set more off
clear


*------------------
* Open cleaned and operationalized survey data
*------------------
use "SL_survey_analysis.dta", clear


	
*------------------	
* Set the data as survey data
*------------------
svyset sl_id
	

		
		

	
	
*-------------------------------------------------------------------------------
* D6 (ordered): have inter-ethnic relationships improved since 2015?
*				1 = much worse 		
*				2 = worse
*				3 = same
*				4 = better
*				5 = much better
*-------------------------------------------------------------------------------


		*-------------------
		* Bivariate Model
		*-------------------
		#delimit ;
		svy:   	ologit
				d6_ethnic_rel_since_2015			
				a3_ethnic_sinhalese
				a3_ethnic_tamil_moors 
		;
		#delimit cr
		estimate store model_d6_biv


		
		*-------------------
		* Full Model
		*-------------------		
		#delimit ;
		svy:   	ologit
				d6_ethnic_rel_since_2015			
				a3_ethnic_sinhalese
				a3_ethnic_tamil_moors 
				x4_gender 				
				a1_age
				age2
				education_low 
				education_high				
				a7_econ_own_condition
				a9_econ_own_vs_all
				e4_phys_harm_war
				c1_c2_consumption_traditional 	// Baseline: family 
				c1_c2_consumption_online
		;
		#delimit cr
		estimate store model_d6



		*-------------------
		* Predicted margins: Sinhalese = 1, all other variables to the mean 
		*-------------------
		
			* Outcome 1 = much worse
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 1	
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 ) 		
					  predict(outcome(1))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_sinhalese_o1
			parmest, saving("ethnicity2015_sinhalese_1.dta", replace) level(95) idstr("sinhalese_1")

			*
			quietly estimates restore model_d6
			
			
			* Outcome 2 = worse
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 1
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 )		
					  predict(outcome(2))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_sinhalese_o2
			parmest, saving("ethnicity2015_sinhalese_2.dta", replace) level(95) idstr("sinhalese_2")
			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 3 = same
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 1
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 )		
					  predict(outcome(3))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_sinhalese_o3
			parmest, saving("ethnicity2015_sinhalese_3.dta", replace) level(95) idstr("sinhalese_3")

			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 4 = better
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 1
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 )		
					  predict(outcome(4))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_sinhalese_o4
			parmest, saving("ethnicity2015_sinhalese_4.dta", replace) level(95) idstr("sinhalese_4")
			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 5 = much better
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 1
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 )		
					  predict(outcome(5))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_sinhalese_o5
			parmest, saving("ethnicity2015_sinhalese_5.dta", replace) level(95) idstr("sinhalese_5")
			*
			quietly estimates restore model_d6
			
		

		*-------------------
		* Predicted margins: Tamil (all ehtnics to = 0), all other variables to the mean 
		*-------------------
		
			* Outcome 1 = much worse
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0	
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 ) 		
					  predict(outcome(1))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tamil_o1
			parmest, saving("ethnicity2015_tamil_1.dta", replace) level(95) idstr("tamil_1")
			*
			quietly estimates restore model_d6
			
			
			* Outcome 2 = worse
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 )		
					  predict(outcome(2))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tamil_o2
			parmest, saving("ethnicity2015_tamil_2.dta", replace) level(95) idstr("tamil_2")
			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 3 = same
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 )		
					  predict(outcome(3))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tamil_o3
			parmest, saving("ethnicity2015_tamil_3.dta", replace) level(95) idstr("tamil_3")
			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 4 = better
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 )		
					  predict(outcome(4))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tamil_o4
			parmest, saving("ethnicity2015_tamil_4.dta", replace) level(95) idstr("tamil_4")
			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 5 = much better
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 )		
					  predict(outcome(5))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tamil_o5
			parmest, saving("ethnicity2015_tamil_5.dta", replace) level(95) idstr("tamil_5")
			*
			quietly estimates restore model_d6


		
		
		*-------------------
		* Predicted margins: Tamil Moors = 1, all other variables to the mean 
		*-------------------
		
			* Outcome 1 = much worse
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0	
						 a3_ethnic_tamil_moors = 1
						 (means) _all
						 ) 		
					  predict(outcome(1))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tm_o1
			parmest, saving("ethnicity2015_tamil_moor_1.dta", replace) level(95) idstr("tamil_moor_1")
			*
			quietly estimates restore model_d6
			
			
			* Outcome 2 = worse
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 1
						 (means) _all
						 )		
					  predict(outcome(2))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tm_o2
			parmest, saving("ethnicity2015_tamil_moor_2.dta", replace) level(95) idstr("tamil_moor_2")
			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 3 = same
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 1
						 (means) _all
						 )		
					  predict(outcome(3))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tm_o3
			parmest, saving("ethnicity2015_tamil_moor_3.dta", replace) level(95) idstr("tamil_moor_3")
			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 4 = better
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 1
						 (means) _all
						 )		
					  predict(outcome(4))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tm_o4
			parmest, saving("ethnicity2015_tamil_moor_4.dta", replace) level(95) idstr("tamil_moor_4")
			*
			quietly estimates restore model_d6
			
			
			
			* Outcome 5 = much better
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 1
						 (means) _all
						 )		
					  predict(outcome(5))  
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d6_tm_o5
			parmest, saving("ethnicity2015_tamil_moor_5.dta", replace) level(95) idstr("tamil_moor_5")
			*
			quietly estimates restore model_d6
				
		
		*------------------
		* Table the predicted probabilities: point estimates, CI in squared brackets
		*------------------		
		#delimit ;
		estout model_d6_sinhalese_o* model_d6_tamil_o* model_d6_tm_o*, cells(b(star fmt(%9.3f))   ci(par([  ,  ]) fmt(2)))  
								 starlevels(+ 0.10 * 0.05 ** 0.01)
								 legend 
								 label collabels(none) 
								 varlabels(_cons "Predicted Probability")
		;
		#delimit cr


		#delimit ;
		estout model_d6_tm_o*, cells(b(star fmt(%9.3f))   ci(par([  ,  ]) fmt(2)))  
								 starlevels(+ 0.10 * 0.05 ** 0.01)
								 legend 
								 label collabels(none) 
								 varlabels(_cons "Predicted Probability")
		;
		#delimit cr

		
		
		
		*------------------
		* Plot the predicted probabilities
		*------------------
		*Preserving data
               
			   preserve
               
               *Appending exported datasets with marginal effects 
               clear
               use "ethnicity2015_sinhalese_1.dta"
               append using "ethnicity2015_sinhalese_2.dta"
			   append using "ethnicity2015_sinhalese_3.dta"
			   append using "ethnicity2015_sinhalese_4.dta"
			   append using "ethnicity2015_sinhalese_5.dta"
			   append using "ethnicity2015_tamil_1.dta"
			   append using "ethnicity2015_tamil_2.dta"
			   append using "ethnicity2015_tamil_3.dta"
			   append using "ethnicity2015_tamil_4.dta"
			   append using "ethnicity2015_tamil_5.dta"

			   
              *Reshaping data
			   gen event = substr(idstr,-1,.)
			   destring event, replace
			   gen id = "Tamil" if regexm(idstr, "tamil")
			   replace id = "Sinhalese" if regexm(idstr, "sinhalese")
			   replace id = "Tamil Moor" if regexm(idstr, "tamil_moor")
               keep estimate min* max* id event 

               
               *Figure: Ethnicity           
               #delimit ;                           
               twoway rcap min95 max95 event if id == "Sinhalese", lcolor(gs0) lwidth(*1.15)
                                || scatter estimate event if id == "Sinhalese", msymbol(o) msize(*.75) mfcolor(white) mlcolor(gs0) mlwidth(*.8) mlalign(outside)                                                                        
                                || line estimate event if id == "Sinhalese", lcolor(gs0) lwidth(*.9) 
								|| rcap min95 max95 event if id == "Tamil", lcolor(ebblue) lwidth(*1.15)                                                       
                                || scatter estimate event if id == "Tamil", msymbol(D) mcolor(ebblue) msize(*.75) lcolor(black)
								|| line estimate event if id == "Tamil", lcolor(ebblue) lwidth(*.9)  lp(shortdash) 
					 xsize(1.2) ysize(1) 
                     xscale(range(.75 3.25) extend titlegap(+2))
                     yscale(range() extend titlegap(+2))
                     ylabel(0(0.1).5 , valuelabel nogrid angle(0) labsize(*1.6))
					 ymtick(##5)
                     xlabel(1 `""Much" "Worse""'
							2 "Worse"
							3 "Same"
							4 "Better"
							5 `""Much" "Better""', angle(horizontal) nogrid labsize(*1.6))
                     xtitle("", size(*.75))
                     ytitle("Pr(Change in ethnic relations)", size(*1.8))
                     legend(order( 2 "Sinhalese" 
								   4 "Tamil")                                                                         
                                   col(1)
                                   size(*1.6)                                          
                                   keygap(*1.2) 
                                   symxsize(*0.35) 
                                   symysize(*0.5)
								   bm(zero)
                                   region(lcolor(gs0)) bm(l=12 r=3)
                                   position(2) ring(0) 
                                   )
					plotregion(lcolor(black)) 
                    graphregion(color(white) lcolor(white) lalign(outside))
					saving(ethnicity2015, replace)
                ; 
               #delimit cr                         
                      
               
               *Restoring data
               restore

		
		

		
		
		
		
		

*-------------------------------------------------------------------------------
* D10 (binary): should gov cracdown for peace?
*				0 = no
*				1 = yes
*-------------------------------------------------------------------------------


	*-------------
	* Bivariate Model
	*-------------
		#delimit ;
		svy:   	logit
				d10_crackdown_peace			
				a3_ethnic_sinhalese
				a3_ethnic_tamil_moors 
		;
		#delimit cr
		estimate store model_d10_biv


	
	*-------------
	* Full Model
	*-------------
		#delimit ;
		svy:   	logit
				d10_crackdown_peace			
				a3_ethnic_sinhalese
				a3_ethnic_tamil_moors 
				x4_gender 				
				a1_age
				age2
				education_low 
				education_high				
				a7_econ_own_condition
				a9_econ_own_vs_all
				e4_phys_harm_war
				c1_c2_consumption_traditional 	// Baseline: family 
				c1_c2_consumption_online
		;
		#delimit cr
		estimate store model_d10
		
		
		

		*-------------------
		* Predicted margins: Sinhalese, tamils, and tamil moors, all other variables to the mean 
		*-------------------
		
			* Outcome 1 = Sinhalese
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 1
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 ) 		
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d10_sinhalese_o1
			*
			quietly estimates restore model_d10
			
			
			* Outcome 1 = Tamil 
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 0
						 (means) _all
						 ) 		
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d10_tamil_o1
			*
			quietly estimates restore model_d10
			
			
			* Outcome 1 = Tamil Moor 
			#delimit ;
			margins,  at(a3_ethnic_sinhalese = 0
						 a3_ethnic_tamil_moors = 1
						 (means) _all
						 ) 		
					  contrast
					  post
			;
			#delimit cr
			estimates store model_d10_tm_o1
			*
			quietly estimates restore model_d10
			
		
		
		*------------------
		* Table the predicted probabilities: point estimates, CI in squared brackets
		*------------------		
		#delimit ;
		estout model_d10_sinhalese_o* model_d10_tamil_o* model_d10_tm_o1, cells(b(star fmt(%9.3f))   ci(par([  ,  ]) fmt(2)))  
								 starlevels(+ 0.10 * 0.05 ** 0.01)
								 legend 
								 label collabels(none) 
								 varlabels(_cons "Predicted Probability")
		;
		#delimit cr

		
		*------------------
		* Plot the predicted probabilities
		*------------------
		
		#delimit ;
		coefplot (model_d10_sinhalese_o1, 	label()  msymbol(O) mlcolor(gs0) 		msize(1.2) mlwidth(.2) mfcolor(gs16) 		ciopts(recast(rcap) lwidth(*1) lcolor(gs0)))
		         (model_d10_tamil_o1, 		label()  msymbol(D) mlcolor(ebblue) 	msize(1.2) mlwidth(.2) mfcolor(ebblue) 	ciopts(recast(rcap) lwidth(*1) lcolor(ebblue)))
				,bycoefs 
				 vertical
				 levels(95)
				 xscale(range(.65 1.3) extend titlegap(+2))
				 ytitle("Pr(Government crackdown)", size(*1.8))
				 ylabel(0(0.1).5, labsize(*1.6) angle(horizontal))
				 ymtick(##5)
				 xtitle("", linegap(2) margin(medium) size(small) height(5))
				 xlabel(0.85 `""Sinhalese" "   ""'  1.2 `""Tamil" "' , labsize(*1.6))
				 ysize (1) xsize(1.2) 
				 scheme(s1mono)
				 legend(off)
				 graphregion(c(white) margin(small))
				 saving(crackdown, replace)
		; 
		#delimit cr	

	   
	   
	   
	  	
	

*-------------------------------------------------------------------------------
* Combine plots: d6 and d10 (Figure 3)
*-------------------------------------------------------------------------------

# delimit ;
gr combine ethnicity2015.gph crackdown.gph, 
	cols(2)
	title("", 
	       ring(1) position(9) orientation(vertical) size(small) linegap(2)
		   margin(r=2))
	subtitle("", 
	       ring(1) position(6) orientation(horizontal) size(small) linegap(2)
		   margin(r=2))
	ysize(1) xsize(2.5)
	plotregion(margin(small) lc(black) ls(solid) lw(vthin) col(white)) 
	graphregion(col(white) lcolor(white) lalign(outside))
	scheme(s2mono)
	imargin(0 8 0 0)
;
#delimit cr

graph display, ysize(1) xsize(2.25) 
graph export "figure_3.pdf", as(pdf) replace		


erase ethnicity2015.gph 
erase crackdown.gph

	
*-------------------------------------------------------------------------------
* Export Appendix table with bivariate + full models for peace, security, speech, and trust only)
*-------------------------------------------------------------------------------

#delimit;													
	esttab 	model_d6
			model_d10
			using "table_3.tex", replace 				
	b(3) se(3) 
	stats(	N,  fmt(0) 
			labels("Number of observations")
			) 
	order( 		a3_ethnic_sinhalese			
				a3_ethnic_tamil_moors 
				x4_gender 					
				a1_age
				age2
				education_low 
				education_high				
				a7_econ_own_condition
				a9_econ_own_vs_all
				e4_phys_harm_war
				c1_c2_consumption_traditional  
				c1_c2_consumption_online
			) 
	coeflabels( a3_ethnic_sinhalese				"Sinhalese"			
				a3_ethnic_tamil_moors 			"Tamil-Moor"
				x4_gender 						"Female"
				a1_age							"Age"
				age2							"Age\textsuperscript{2}"
				education_low 					"Low education"
				education_high					"High education"
				a7_econ_own_condition			"Absolute economic situation"
				a9_econ_own_vs_all				"Relative economic situation"
				e4_phys_harm_war				"Harmed in war"
				c1_c2_consumption_traditional  	"Use traditional media"
				c1_c2_consumption_online		"Use online media"
			) 
	label 
	star(\dag 0.10 * 0.05  ** 0.01  *** 0.001) 
	nonumber
	mgroups("\multicolumn{1}{c}{(6)}" 
			"\multicolumn{1}{c}{(7)}" 
			, pattern(1 1))
	nobaselevels 
	eqlabels(none)
	nomtitles
	nolegend 
	nonotes 
	booktabs 
	alignment(D{.}{.}{-1}) 
	compress
	substitute(\_ _)
	fragment
	prehead(
	\begin{table}[!ht]\centering
	\def\sym#1{\ifmmode^{#1}\else\(^{#1}\)\fi}
		\caption{Perceived changes in ethnic relations since 2015 and preference for government harsh responses \label{tab:table3}}
		\setlength{\tabcolsep}{2pt}
		\scalebox{.75}{		
		\begin{tabular}{l*{2}{S}}
		\hline	
		\addlinespace
		& \multicolumn{2}{c}{\textit{Dependent variable}}\\ 
		\cmidrule(r{5pt}l{5pt}){2-3} 
        \addlinespace
		& \multicolumn{1}{c}{Have inter-ethnic relations } & \multicolumn{1}{c}{Government crackdown on } \\ 
		& \multicolumn{1}{c}{improved since 2015} & \multicolumn{1}{c}{extremists would promote peace} \\
		)
		postfoot(
		\hline
		\multicolumn{3}{l}{\footnotesize \textit{Note:} Coefficients with robust standard errors in parentheses.}\\
		\multicolumn{3}{l}{\footnotesize  $\dag$ p$<$0.1, * p$<$0.05, ** p$<$0.01, *** p$<$0.001.}\\
		\end{tabular}
		}
		\end{table}		
		);
#delimit cr	
	
	
	
	

	
		
*-------------------------------------------------------------------------------
* Export Appendix table with bivariate + full models for peace, security, speech, and trust only)
*-------------------------------------------------------------------------------

#delimit;													
	esttab 	model_d6_biv
			model_d6
			model_d10_biv
			model_d10
			using "table_a5.tex", replace 				
	b(3) se(3) 
	stats(	N,  fmt(0) 
			labels("Number of observations")
			) 
	order( 		a3_ethnic_sinhalese			
				a3_ethnic_tamil_moors 
				x4_gender 					
				a1_age
				age2
				education_low 
				education_high				
				a7_econ_own_condition
				a9_econ_own_vs_all
				e4_phys_harm_war
				c1_c2_consumption_traditional  
				c1_c2_consumption_online
			) 
	coeflabels( a3_ethnic_sinhalese				"Sinhalese"			
				a3_ethnic_tamil_moors 			"Tamil-Moor"
				x4_gender 						"Female"
				a1_age							"Age"
				age2							"Age\textsuperscript{2}"
				education_low 					"Low education"
				education_high					"High education"
				a7_econ_own_condition			"Absolute economic situation"
				a9_econ_own_vs_all				"Relative economic situation"
				e4_phys_harm_war				"Harmed in war"
				c1_c2_consumption_traditional  	"Use traditional media"
				c1_c2_consumption_online		"Use online media"
			) 
	label 
	star(\dag 0.10 * 0.05  ** 0.01  *** 0.001) 
	nonumber
	mgroups("\multicolumn{1}{c}{(1)}" 
			"\multicolumn{1}{c}{(2)}" 
			"\multicolumn{1}{c}{(3)}"
			"\multicolumn{1}{c}{(4)}"
			, pattern(1 1 1 1 ))
	nobaselevels 
	eqlabels(none)
	nomtitles
	nolegend 
	nonotes 
	booktabs 
	alignment(D{.}{.}{-1}) 
	compress
	substitute(\_ _)
	fragment
	prehead(
	\begin{table}[!ht]\centering
	\def\sym#1{\ifmmode^{#1}\else\(^{#1}\)\fi}
		\caption{Perceived changes in ethnic relations since 2015 and preference for government harsh responses \label{tab:table_rob_67}}
		\setlength{\tabcolsep}{2pt}
		\scalebox{.75}{		
		\begin{tabular}{l*{4}{S}}
		\hline	
		\addlinespace
		& \multicolumn{4}{c}{\textit{Dependent variable}}\\ 
		\cmidrule(r{5pt}l{5pt}){2-5} 
        \addlinespace
		& \multicolumn{2}{c}{Have inter-ethnic relations } & \multicolumn{2}{c}{Government crackdown on } \\ 
		& \multicolumn{2}{c}{improved since 2015} & \multicolumn{2}{c}{extremists would promote peace} \\ 
		\cmidrule(r{5pt}l{5pt}){2-3} \cmidrule(r{5pt}l{5pt}){4-5} 
		)
		postfoot(
		\hline
		\multicolumn{5}{l}{\footnotesize \textit{Note:} Coefficients with robust standard errors in parentheses.}\\
		\multicolumn{5}{l}{\footnotesize  $\dag$ p$<$0.1, * p$<$0.05, ** p$<$0.01, *** p$<$0.001.}\\
		\end{tabular}
		}
		\end{table}		
		);
#delimit cr


